# Check requisite packages are installed.
packages <- c(
"plotly",
"dplyr"
)
for (pkg in packages) {
library(pkg, character.only = TRUE)
}
package 㤼㸱plotly㤼㸲 was built under R version 4.0.5Loading required package: ggplot2
package 㤼㸱ggplot2㤼㸲 was built under R version 4.0.3
Attaching package: 㤼㸱plotly㤼㸲
The following object is masked from 㤼㸱package:ggplot2㤼㸲:
last_plot
The following object is masked from 㤼㸱package:stats㤼㸲:
filter
The following object is masked from 㤼㸱package:graphics㤼㸲:
layout
package 㤼㸱dplyr㤼㸲 was built under R version 4.0.4
Attaching package: 㤼㸱dplyr㤼㸲
The following objects are masked from 㤼㸱package:stats㤼㸲:
filter, lag
The following objects are masked from 㤼㸱package:base㤼㸲:
intersect, setdiff, setequal, union
Pulling code almost directly from LM1996-NumPoolComScaling-Results-2021-05.Rmd.
dirViking <- c(
file.path(
getwd(), "LCAB_LawMorton1996-NumericalPoolCommunityScaling"
),
file.path(
getwd(), "LCAB_LawMorton1996-NumericalPoolCommunityScaling2"
)
)
dirVikingResults <- file.path(
dirViking, c("results-2021-04", "save-2021-05-10") # Latter not 100% yet.
)
resultFormat <- paste0(
"run-",
"%d", # Combination Number, or CombnNum.
"-",
"%s", # Run Seed.
".RDS"
)
# Copied from LawMorton1996-NumericalPoolCommunityScaling-Calculation.R
set.seed(38427042)
basal <- c(3, 10, 30, 100, 300, 1000)
consumer <- c(3, 10, 30, 100, 300, 1000) * 2
events <- (max(basal) + max(consumer)) * 2
runs <- 100
logBodySize <- c(-2, -1, -1, 1) # Morton and Law 1997 version.
parameters <- c(0.01, 10, 0.5, 0.2, 100, 0.1)
# Need to rerun seedsPrep to get the random number generation right for seedsRun
seedsPrep <- runif(2 * length(basal) * length(consumer)) * 1E8
seedsRun <- runif(runs * length(basal) * length(consumer)) * 1E8
paramFrame <- with(list(
b = rep(basal, times = length(consumer)),
c = rep(consumer, each = length(basal)),
s1 = seedsPrep[1:(length(basal) * length(consumer))],
s2 = seedsPrep[
(length(basal) * length(consumer) + 1):(
2 * length(basal) * length(consumer))
],
sR = seedsRun
), {
temp <- data.frame(
CombnNum = 0,
Basals = b,
Consumers = c,
SeedPool = s1,
SeedMat = s2,
SeedRuns = "",
SeedRunsNum = 0,
EndStates = I(rep(list(""), length(b))),
EndStatesNum = 0,
EndStateSizes = I(rep(list(""), length(b))),
EndStateSizesNum = NA,
EndStateAssembly = I(rep(list(""), length(b))),
EndStateAbundance = I(rep(list(""), length(b))),
Dataset = "2021-04",
DatasetID = 1,
stringsAsFactors = FALSE
)
for (i in 1:nrow(temp)) {
seeds <- sR[((i - 1) * runs + 1) : (i * runs)]
temp$SeedRuns[i] <- toString(seeds) # CSV
temp$SeedRunsNum[i] <- length(seeds)
}
temp$CombnNum <- 1:nrow(temp)
temp
})
# Note: n + 2 end states. Failure to finish, failure to obtain state, and state.
for (i in 1:nrow(paramFrame)) {
resultsList <- list(
"No Run" = 0,
"No State" = 0
)
resultsSize <- list(
"0" = 0
)
resultsAssembly <- list(
"No Run" = data.frame(),
"No State" = data.frame()
)
seeds <- unlist(strsplit(paramFrame$SeedRuns[i], ', '))
for (seed in seeds) {
fileName <- file.path(
dirVikingResults[paramFrame$DatasetID[i]],
sprintf(resultFormat, paramFrame$CombnNum[i], seed)
)
if (file.exists(fileName)) {
temp <- load(fileName)
temp <- eval(parse(text = temp)) # Get objects.
if (is.data.frame(temp)) {
community <- toString(
temp[[ncol(temp)]][[nrow(temp)]]
)
size <- toString(length(temp[[ncol(temp)]][[nrow(temp)]]))
if (community == "") {
resultsList$`No State` <- resultsList$`No State` + 1
resultsSize$`0` <- resultsSize$`0` + 1
} else if (community %in% names(resultsList)) {
resultsList[[community]] <- resultsList[[community]] + 1
resultsSize[[size]] <- resultsSize[[size]] + 1
} else {
resultsList[[community]] <- 1
resultsAssembly[[community]] <- temp
if (size %in% resultsSize) {
resultsSize[[size]] <- resultsSize[[size]] + 1
} else {
resultsSize[[size]] <- 1
}
}
} else {
resultsList$`No State` <- resultsList$`No State` + 1
resultsSize$`0` <- resultsSize$`0` + 1
}
} else {
resultsList$`No Run` <- resultsList$`No Run` + 1
resultsSize$`0` <- resultsSize$`0` + 1
}
}
paramFrame$EndStates[[i]] <- resultsList
paramFrame$EndStatesNum[i] <- length(resultsList) - 2 # ! No State, No Run
paramFrame$EndStateSizes[[i]] <- resultsSize
paramFrame$EndStateSizesNum[i] <- length(resultsSize) - 1 # ! 0
paramFrame$EndStateAssembly[[i]] <- resultsAssembly
}
source(
file.path(getwd(),
"LawMorton1996-NumericalPoolCommunityScaling-Settings2.R")
)
oldNrow <- nrow(paramFrame)
paramFrame <- rbind(paramFrame, with(list(
b = rep(basal, times = length(consumer)),
c = rep(consumer, each = length(basal)),
s1 = seedsPrep[1:(length(basal) * length(consumer))],
s2 = seedsPrep[
(length(basal) * length(consumer) + 1):(
2 * length(basal) * length(consumer))
],
sR = seedsRun
), {
temp <- data.frame(
CombnNum = 0,
Basals = b,
Consumers = c,
SeedPool = s1,
SeedMat = s2,
SeedRuns = "",
SeedRunsNum = 0,
EndStates = I(rep(list(""), length(b))),
EndStatesNum = 0,
EndStateSizes = I(rep(list(""), length(b))),
EndStateSizesNum = NA,
EndStateAssembly = I(rep(list(""), length(b))),
EndStateAbundance = I(rep(list(""), length(b))),
Dataset = "2021-05",
DatasetID = 2,
stringsAsFactors = FALSE
)
for (i in 1:nrow(temp)) {
seeds <- sR[((i - 1) * runs + 1) : (i * runs)]
temp$SeedRuns[i] <- toString(seeds) # CSV
temp$SeedRunsNum[i] <- length(seeds)
}
temp$CombnNum <- 1:nrow(temp)
temp
})
)
# Note: n + 2 end states. Failure to finish, failure to obtain state, and state.
# Modified from above, but with the abundance recorded.
for (i in (oldNrow + 1):nrow(paramFrame)) {
resultsList <- list(
"No Run" = 0,
"No State" = 0
)
resultsSize <- list(
"0" = 0
)
resultsAssembly <- list(
"No Run" = data.frame(),
"No State" = data.frame()
)
resultsAbund <- list(
"No Run" = "",
"No State" = ""
)
seeds <- unlist(strsplit(paramFrame$SeedRuns[i], ', '))
for (seed in seeds) {
fileName <- file.path(
dirVikingResults[paramFrame$DatasetID[i]],
sprintf(resultFormat, paramFrame$CombnNum[i], seed)
)
if (file.exists(fileName)) {
temp <- load(fileName)
temp <- eval(parse(text = temp)) # Get objects.
if (is.list(temp) && "Result" %in% names(temp)) {
if (is.data.frame(temp$Result))
community <- temp$Result$Community[[nrow(temp$Result)]]
else
community <- temp$Result
size <- toString(length(community))
if (community[1] != "")
abund <- toString(temp$Abund[community + 1])
else
abund <- ""
community <- toString(community)
if (community == "") {
resultsList$`No State` <- resultsList$`No State` + 1
resultsSize$`0` <- resultsSize$`0` + 1
} else if (community %in% names(resultsList)) {
resultsList[[community]] <- resultsList[[community]] + 1
resultsSize[[size]] <- resultsSize[[size]] + 1
} else {
resultsList[[community]] <- 1
resultsAssembly[[community]] <- temp
resultsAbund[[community]] <- abund
if (size %in% resultsSize) {
resultsSize[[size]] <- resultsSize[[size]] + 1
} else {
resultsSize[[size]] <- 1
}
}
} else {
resultsList$`No State` <- resultsList$`No State` + 1
resultsSize$`0` <- resultsSize$`0` + 1
}
} else {
resultsList$`No Run` <- resultsList$`No Run` + 1
resultsSize$`0` <- resultsSize$`0` + 1
}
}
paramFrame$EndStates[[i]] <- resultsList
paramFrame$EndStatesNum[i] <- length(resultsList) - 2 # ! No State, No Run
paramFrame$EndStateSizes[[i]] <- resultsSize
paramFrame$EndStateSizesNum[i] <- length(resultsSize) - 1 # ! 0
paramFrame$EndStateAssembly[[i]] <- resultsAssembly
paramFrame$EndStateAbundance[[i]] <- resultsAbund
}
testRowNums <- nrow(paramFrame)
testRowsToAdd <- c(2, 6) # Make sure to put in numerical order!
paramFrame <- with(
list(
basal2 = c(5, 10, 15),
consumer2 = c(20, 40, 60),
logBodySize = c(-2, -1, -1, 0),
parameters = c(0.01, 10, 0.5, 0.2, 100, 0.1)
),
{
set.seed(3680180)
seedsPrep2 <- runif(2 * length(basal2) * length(consumer2)) * 1E8
with(list(
b = rep(basal2, times = length(consumer2)),
c = rep(consumer2, each = length(basal2)),
s1 = seedsPrep2[1:(length(basal2) * length(consumer2))],
s2 = seedsPrep2[
(length(basal2) * length(consumer2) + 1):(
2 * length(basal2) * length(consumer2))
]
), {
rbind(
paramFrame,
data.frame(
CombnNum = testRowsToAdd,
Basals = b[testRowsToAdd],
Consumers = c[testRowsToAdd],
SeedPool = s1[testRowsToAdd],
SeedMat = s2[testRowsToAdd],
SeedRuns = "",
SeedRunsNum = 0,
EndStates = I(rep(list(""), length(testRowsToAdd))),
EndStatesNum = 0,
EndStateSizes = I(rep(list(""), length(testRowsToAdd))),
EndStateSizesNum = NA,
EndStateAssembly = I(rep(list(""), length(testRowsToAdd))),
EndStateAbundance = I(rep(list(""), length(testRowsToAdd))),
Dataset = "Test",
DatasetID = max(paramFrame$DatasetID) + 1,
stringsAsFactors = FALSE
)
)
}
)
}
)
testRowNums <- (testRowNums + 1):nrow(paramFrame)
resultsList <- list(
list(
"No Run" = 0,
"No State" = 0,
"2, 4, 6, 12, 29" = 1,
"2, 4, 6, 13, 29" = 1
),
list(
"No Run" = 0,
"No State" = 0,
"8, 10, 12, 14, 15, 16, 39, 43" = 1,
"8, 12, 14, 15, 16, 38, 39" = 1
)
)
resultsSize <- list(
list(
"0" = 0,
"5" = 2
),
list(
"0" = 0,
"8" = 1,
"7" = 1
)
)
resultsAbund <- list(
list(
"No Run" = "",
"No State" = "",
"2, 4, 6, 12, 29" = "742.88553671712, 80.579233072626, 162.128399850253, 20.2082198699389, 18.8589490510429",
"2, 4, 6, 13, 29" = "668.664143581837, 119.024146851052, 127.680269383867, 30.657960866033, 13.4844194707944"
),
list(
"No Run" = "",
"No State" = "",
"8, 10, 12, 14, 15, 16, 39, 43" = "20.7665807606491, 32.4461165261454, 80.4033387818895, 817.879722033354, 121.136570782828, 18.0390671088957, 12.3834561177271, 19.9674718543196",
"8, 12, 14, 15, 16, 38, 39" = "82.592048492812, 138.267379166014, 938.158436379166, 51.8610963745021, 5.03556251837491, 14.1019343145825, 25.9231062711228"
)
)
for (j in seq_along(testRowNums)) {
i <- testRowNums[j]
paramFrame$EndStates[[i]] <- resultsList[[j]]
paramFrame$EndStatesNum[i] <- length(resultsList[[j]]) - 2 # ! No State, No Run
paramFrame$EndStateSizes[[i]] <- resultsSize[[j]]
paramFrame$EndStateSizesNum[i] <- length(resultsSize[[j]]) - 1 # ! 0
paramFrame$EndStateAbundance[[i]] <- resultsAbund[[j]]
}
# X, Y, Basal and Consumer.
# Z = Sizes of the Endstates.
plotScalingData <- data.frame(
CombnNum = rep(paramFrame$CombnNum, paramFrame$EndStatesNum),
Basals = rep(paramFrame$Basals, paramFrame$EndStatesNum),
Consumers = rep(paramFrame$Consumers, paramFrame$EndStatesNum),
Dataset = rep(paramFrame$Dataset, paramFrame$EndStatesNum),
DatasetID = rep(paramFrame$DatasetID, paramFrame$EndStatesNum)
)
# Communities
comms <- unlist(lapply(paramFrame$EndStates, names))
freqs <- unlist(paramFrame$EndStates)
asmbl <- unlist(paramFrame$EndStateAssembly, recursive = FALSE)
asmbl <- asmbl[comms != "No Run" & comms != "No State"]
freqs <- freqs[comms != "No Run" & comms != "No State"]
comms <- comms[comms != "No Run" & comms != "No State"]
asmbl <- lapply(asmbl, function(d) {
if (is.null(d)) return(NA)
if ("Result.Outcome" %in% names(d))
d %>% dplyr::filter(Result.Outcome != "Type 1 (Failure)" &
Result.Outcome != "Present")
else
d$Result %>% dplyr::filter(Outcome != "Type 1 (Failure)" &
Outcome != "Present")
})
plotScalingData$Communities <- comms
plotScalingData$CommunityFreq <- freqs
plotScalingData$CommunitySeq <- asmbl
# Community Size
temp <- unlist(lapply(strsplit(plotScalingData$Communities, ','), length))
plotScalingData$CommunitySize <- temp
# For usage by the reader.
plotScaling <- plotly::plot_ly(
plotScalingData,
x = ~Basals,
y = ~Consumers,
z = ~CommunitySize,
color = ~Dataset,
colors = c("red", "blue", "black")
)
plotScaling <- plotly::add_markers(plotScaling)
plotScaling <- plotly::layout(
plotScaling,
scene = list(
xaxis = list(type = "log"),
yaxis = list(type = "log"),
camera = list(
eye = list(
x = -1.25, y = -1.25, z = .05
)
)
)
)
plotScaling
# > runif(1) * 1E8
# [1] 82598679
set.seed(82598679)
mats <- list()
poolsall <- list() # name pools used in save data; be careful!
for (i in 1:length(dirViking)) {
temp <- load(file.path(
dirViking[i],
paste0("LawMorton1996-NumericalPoolCommunityScaling-PoolMats",
if (i > 1) i else "",
".RDS")
))
mats[[i]] <- eval(parse(text = temp[1]))
poolsall[[i]] <- eval(parse(text = temp[2]))
}
pools <- poolsall
# Add in the test datasets.
poolsTemp <- list()
matsTemp <- list()
for (r in testRowNums) {
testRowRow <- paramFrame[r, ]
poolsTemp[[testRowRow$CombnNum]] <- with(testRowRow,
RMTRCode2::LawMorton1996_species(
Basal = Basals,
Consumer = Consumers,
Parameters = c(0.01, 10, 0.5, 0.2, 100, 0.1),
LogBodySize = c(-2, -1, -1, 0),
seed = SeedPool
)
)
matsTemp[[testRowRow$CombnNum]] <- with(testRowRow,
RMTRCode2::LawMorton1996_CommunityMat(
Pool = poolsTemp[[CombnNum]],
Parameters = c(0.01, 10, 0.5, 0.2, 100, 0.1),
seed = SeedMat
)
)
}
executing %dopar% sequentially: no parallel backend registered
pools[[i + 1]] <- poolsTemp
mats[[i + 1]] <- matsTemp
oldCandidateData <- load(file.path(getwd(), "candidateDataSoFar.Rdata"))
oldCandidateData <- eval(parse(text = oldCandidateData))
candidateData <- plotScalingData %>% dplyr::group_by(
CombnNum, Dataset
) %>% dplyr::mutate(
OtherSteadyStates = dplyr::n() - 1
) %>% dplyr::filter(
OtherSteadyStates > 0
)
candidateData %>% dplyr::select(-CommunitySeq)
# First, check if it is in the paramFrame.
# Second, check if it is in the saved data from the previous.
# Otherwise, ignore it, we'll figure out what it is and why it is missing later.
candidateData$CommunityAbund <- ""
for (r in 1:nrow(candidateData)) {
# ID 1:4 are used to identify paramFrame, 5 used to identify abundance
ID <- candidateData[r, 1:6]
paramFrameRow <- paramFrame %>% dplyr::filter(
CombnNum == ID$CombnNum,
Basals == ID$Basals,
Consumers == ID$Consumers,
Dataset == ID$Dataset
)
if (is.list(paramFrameRow$EndStateAbundance[[1]])) {
entry <- which(ID$Communities == names(paramFrameRow$EndStateAbundance[[1]]))
if (length(entry)) {
candidateData$CommunityAbund[r] <- paramFrameRow$EndStateAbundance[[1]][[entry]]
next()
}
}
if (ID$Dataset == "2021-04") {
oldCandDatRow <- oldCandidateData %>% dplyr::filter(
CombnNum == ID$CombnNum,
Basals == ID$Basals,
Consumers == ID$Consumers,
Communities == ID$Communities
)
if (nrow(oldCandDatRow) > 0) {
if (oldCandDatRow$CommunityAbund != "") {
candidateData$CommunityAbund[r] <- oldCandDatRow$CommunityAbund
}
}
}
}
for (r in 1:nrow(candidateData)) {
if (!(candidateData$CommunityAbund[r] == "Failure" |
candidateData$CommunityAbund[r] == "")) next
# Random guesses, starting from structured.
temp <- with(
candidateData[r, ],
RMTRCode2::FindSteadyStateFromEstimate(
Pool = pools[[DatasetID]][[CombnNum]],
InteractionMatrix = mats[[DatasetID]][[CombnNum]],
Community = Communities,
Populations = ifelse(
pools[[DatasetID]][[CombnNum]]$Type[
RMTRCode2::CsvRowSplit(Communities)
] == "Basal",
1000, 10)
)
)
if (any(temp) < 1E-4) {
temp <- "EstimateFailure"
} else {
temp <- toString(temp)
}
candidateData$CommunityAbund[r] <- temp
}
DLSODA- At current T (=R1), MXSTEP (=I1) steps
taken on this call before reaching TOUT
In above message, I1 = 5000
In above message, R1 = 2843.06
an excessive amount of work (> maxsteps ) was done, but integration was not successful - increase maxstepsReturning early. Results are accurate, as far as they go
DLSODA- At current T (=R1), MXSTEP (=I1) steps
taken on this call before reaching TOUT
In above message, I1 = 5000
In above message, R1 = 2876.33
an excessive amount of work (> maxsteps ) was done, but integration was not successful - increase maxstepsReturning early. Results are accurate, as far as they go
DLSODA- At current T (=R1), MXSTEP (=I1) steps
taken on this call before reaching TOUT
In above message, I1 = 5000
In above message, R1 = 2870.65
an excessive amount of work (> maxsteps ) was done, but integration was not successful - increase maxstepsReturning early. Results are accurate, as far as they go
DLSODA- At current T (=R1), MXSTEP (=I1) steps
taken on this call before reaching TOUT
In above message, I1 = 5000
In above message, R1 = 2873.23
an excessive amount of work (> maxsteps ) was done, but integration was not successful - increase maxstepsReturning early. Results are accurate, as far as they go
DLSODA- At current T (=R1), MXSTEP (=I1) steps
taken on this call before reaching TOUT
In above message, I1 = 5000
In above message, R1 = 2871.69
an excessive amount of work (> maxsteps ) was done, but integration was not successful - increase maxstepsReturning early. Results are accurate, as far as they gocoercing argument of type 'double' to logicalcoercing argument of type 'double' to logicalcoercing argument of type 'double' to logicalcoercing argument of type 'double' to logicalcoercing argument of type 'double' to logicalcoercing argument of type 'double' to logical
DLSODA- At current T (=R1), MXSTEP (=I1) steps
taken on this call before reaching TOUT
In above message, I1 = 5000
In above message, R1 = 6597.38
an excessive amount of work (> maxsteps ) was done, but integration was not successful - increase maxstepsReturning early. Results are accurate, as far as they gocoercing argument of type 'double' to logical
DLSODA- At current T (=R1), MXSTEP (=I1) steps
taken on this call before reaching TOUT
In above message, I1 = 5000
In above message, R1 = 3592.76
an excessive amount of work (> maxsteps ) was done, but integration was not successful - increase maxstepsReturning early. Results are accurate, as far as they go
DLSODA- At current T (=R1), MXSTEP (=I1) steps
taken on this call before reaching TOUT
In above message, I1 = 5000
In above message, R1 = 3780.53
an excessive amount of work (> maxsteps ) was done, but integration was not successful - increase maxstepsReturning early. Results are accurate, as far as they go
DLSODA- At current T (=R1), MXSTEP (=I1) steps
taken on this call before reaching TOUT
In above message, I1 = 5000
In above message, R1 = 3786.89
an excessive amount of work (> maxsteps ) was done, but integration was not successful - increase maxstepsReturning early. Results are accurate, as far as they go
DLSODA- At current T (=R1), MXSTEP (=I1) steps
taken on this call before reaching TOUT
In above message, I1 = 5000
In above message, R1 = 3784.07
an excessive amount of work (> maxsteps ) was done, but integration was not successful - increase maxstepsReturning early. Results are accurate, as far as they go
DLSODA- At current T (=R1), MXSTEP (=I1) steps
taken on this call before reaching TOUT
In above message, I1 = 5000
In above message, R1 = 3778.43
an excessive amount of work (> maxsteps ) was done, but integration was not successful - increase maxstepsReturning early. Results are accurate, as far as they go
DLSODA- At current T (=R1), MXSTEP (=I1) steps
taken on this call before reaching TOUT
In above message, I1 = 5000
In above message, R1 = 3777.75
an excessive amount of work (> maxsteps ) was done, but integration was not successful - increase maxstepsReturning early. Results are accurate, as far as they go
DLSODA- At current T (=R1), MXSTEP (=I1) steps
taken on this call before reaching TOUT
In above message, I1 = 5000
In above message, R1 = 3774.98
an excessive amount of work (> maxsteps ) was done, but integration was not successful - increase maxstepsReturning early. Results are accurate, as far as they go
DLSODA- At current T (=R1), MXSTEP (=I1) steps
taken on this call before reaching TOUT
In above message, I1 = 5000
In above message, R1 = 3781.84
an excessive amount of work (> maxsteps ) was done, but integration was not successful - increase maxstepsReturning early. Results are accurate, as far as they gocoercing argument of type 'double' to logical
candidateData <- candidateData %>% dplyr::filter(CommunityAbund != "",
CommunityAbund != "Failure",
CommunityAbund != "EstimateFailure")
candidateData$CommunityProd <- NA
for (r in 1:nrow(candidateData)) {
candidateData$CommunityProd[r] <- with(
candidateData[r, ],
RMTRCode2::Productivity(
Pool = pools[[DatasetID]][[CombnNum]],
InteractionMatrix = mats[[DatasetID]][[CombnNum]],
Community = Communities,
Populations = CommunityAbund
)
)
}
islandFUN <- function(i, dat, pool, mat, dmat) {
temp <- dat[i, ]
RMTRCode2::IslandDynamics(
Pool = pool,
InteractionMatrix = mat,
Communities = c(
list(temp$Communities[1]),
rep("", nrow(dmat) - 2),
temp$Communities[2]
),
Populations = c(
list(temp$CommunityAbund[1]),
rep("", nrow(dmat) - 2),
list(temp$CommunityAbund[2])
),
DispersalPool = 0.0001,
DispersalIsland = dmat,
Verbose = FALSE
)
}
# For each group-dataset,
# For each pair,
# Run Island Dynamics,
# Save the result with its pairing
candidateData$TotalID <- paste(candidateData$CombnNum, candidateData$DatasetID)
islandInteractionsOneTwo <- list()
for (grp in unique(candidateData$TotalID)) {
candidateDataSubset <- candidateData %>% dplyr::filter(TotalID == grp)
if (nrow(candidateDataSubset) == 1) next()
pairingResults <- combn(
nrow(candidateDataSubset), 2,
islandFUN,
dat = candidateDataSubset,
pool = pools[[
candidateDataSubset$DatasetID[1]
]][[candidateDataSubset$CombnNum[1]]],
mat = mats[[
candidateDataSubset$DatasetID[1]
]][[candidateDataSubset$CombnNum[1]]],
dmat = matrix(c(0, 1, 1, 0), nrow = 2, ncol = 2),
simplify = FALSE
)
pairingResults <- lapply(
pairingResults, function(mat, isles) {
mat <- mat[nrow(mat), -1]
retVal <- list()
species <- length(mat) / isles
for (i in 1:isles) {
retVal[[i]] <- mat[((i - 1) * species + 1) : (i * species)]
}
retVal
},
isles = 2
)
islandInteractionsOneTwo[[grp]] <- pairingResults
}
islandInteractionsOneEmptyTwo <- list()
for (grp in unique(candidateData$TotalID)) {
candidateDataSubset <- candidateData %>% dplyr::filter(TotalID == grp)
if (nrow(candidateDataSubset) == 1) next()
pairingResults <- combn(
nrow(candidateDataSubset), 2,
islandFUN,
dat = candidateDataSubset,
pool = pools[[
candidateDataSubset$DatasetID[1]
]][[candidateDataSubset$CombnNum[1]]],
mat = mats[[
candidateDataSubset$DatasetID[1]
]][[candidateDataSubset$CombnNum[1]]],
dmat = matrix(c(
0, 1, 0, # Island 2 -> 1
1, 0, 1, # Island 1 -> 2, Island 3 -> 2
0, 1, 0 # Island 2 -> 3
), nrow = 3, ncol = 3, byrow = TRUE),
simplify = FALSE
)
pairingResults <- lapply(
pairingResults, function(mat, isles) {
mat <- mat[nrow(mat), -1]
retVal <- list()
species <- length(mat) / isles
for (i in 1:isles) {
retVal[[i]] <- mat[((i - 1) * species + 1) : (i * species)]
}
retVal
},
isles = 3
)
islandInteractionsOneEmptyTwo[[grp]] <- pairingResults
}
# Format of table should be:
# ID, Community 1, Community 2, Outcomes 1-2, Outcomes 1-0-2
# For outcomes, species presence will be used.
communities <- NULL
totalCommunities <- NULL
for (grp in unique(candidateData$TotalID)) {
candidateDataSubset <- candidateData %>% dplyr::filter(TotalID == grp)
if (nrow(candidateDataSubset) > 1) {
newCommunities <- combn(
candidateDataSubset$Communities, 2,
)
communities <- c(communities, newCommunities)
totalCommunities <- c(
totalCommunities,
toString(sort(unique(unlist(lapply(newCommunities,
RMTRCode2::CsvRowSplit)))))
)
}
}
islandInteractionsOneTwoWhich <- unlist(lapply(
seq_along(islandInteractionsOneTwo), function(i, x, tC) {
lapply(x[[i]], function(y, tC) {
lapply(y, function(z, tC) {
toString(RMTRCode2::CsvRowSplit(tC)[which(z > 1E-6)])
}, tC = tC)
},
tC = tC[i])
},
x = islandInteractionsOneTwo,
tC = totalCommunities))
islandInteractionsOneEmptyTwoWhich <- unlist(lapply(
seq_along(islandInteractionsOneEmptyTwo), function(i, x, tC) {
lapply(x[[i]], function(y, tC) {
lapply(y, function(z, tC) {
toString(RMTRCode2::CsvRowSplit(tC)[which(z > 1E-6)])
}, tC = tC)
},
tC = tC[i])
},
x = islandInteractionsOneEmptyTwo,
tC = totalCommunities))
islandInteractionResults <- data.frame(
DatasetID = rep(names(islandInteractionsOneTwo),
unlist(lapply(islandInteractionsOneTwo, length))),
Community1 = communities[seq(from = 1, to = length(communities), by = 2)],
Community2 = communities[seq(from = 2, to = length(communities), by = 2)],
OutcomeWOEmpty_Island1 = islandInteractionsOneTwoWhich[
seq(from = 1, to = length(islandInteractionsOneTwoWhich), by = 2)],
OutcomeWOEmpty_Island2 = islandInteractionsOneTwoWhich[
seq(from = 1, to = length(islandInteractionsOneTwoWhich), by = 2)],
OutcomeWEmpty_Island1 = islandInteractionsOneEmptyTwoWhich[
seq(from = 1, to = length(islandInteractionsOneEmptyTwoWhich), by = 3)],
OutcomeWEmpty_Island2 = islandInteractionsOneEmptyTwoWhich[
seq(from = 1, to = length(islandInteractionsOneEmptyTwoWhich), by = 3)],
OutcomeWEmpty_Island3 = islandInteractionsOneEmptyTwoWhich[
seq(from = 1, to = length(islandInteractionsOneEmptyTwoWhich), by = 3)]
)
islandInteractionResults
islandInteractionResults %>% dplyr::mutate(
C1WOInvaded = Community1 != OutcomeWOEmpty_Island1,
C2WOInvaded = Community2 != OutcomeWOEmpty_Island2,
C1WInvaded = Community1 != OutcomeWEmpty_Island1,
C2WInvaded = Community2 != OutcomeWEmpty_Island3,
StalemateWO = !C1WOInvaded & !C2WOInvaded,
StalemateW = !C1WInvaded & !C2WInvaded,
HybridWO = C1WOInvaded & C2WOInvaded,
HybridW = C1WInvaded & C2WInvaded,
) %>% dplyr::select(-dplyr::starts_with("Outcome"))
save(
candidateData,
islandInteractionsOneEmptyTwo,
islandInteractionsOneEmptyTwoWhich,
islandInteractionsOneTwo,
islandInteractionsOneTwoWhich,
mats,
paramFrame,
plotScalingData,
pools,
file = "LM1996-NumPoolCom-QDat-2021-05.RData"
)